Delay Induced Instabilities in Self-Propelling Swarms 
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We consider a general model of self-propelling particles interacting through a pairwise attractive 
force in the presence of noise and communication time delay. Previous work by Erdmann, et al. [Phys. 
Rev. E 71, 051904 (2005)] has shown that a large enough noise intensity will cause a translating 
swarm of individuals to transition to a rotating swarm with a stationary center of mass. We show 
that with the addition of a time delay, the model possesses a transition that depends on the size of the 
coupling amplitude. This transition is independent of the initial swarm state (traveling or rotating) 
and is characterized by the alignment of all of the individuals along with a swarm oscillation. By 
considering the mean field equations without noise, we show that the time delay induced transition 
is associated with a Hopf bifurcation. The analytical result yields good agreement with numerical 
computations of the value of the coupling parameter at the Hopf point. 
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The collective motion of multi-agent systems have long 
been observed in biological populations including bacte- 
rial colonies [l], 0, [3] , slime molds 0, [B| , locusts @ and 
fish [?| ■ However, mathematical studies of swarming be- 
havior have been performed for only a few decades. In ad- 
dition to providing examples of biological pattern forma- 
tion, the information gained from these mathematical in- 
vestigations has led to an increased ability to int ellig ently 
design and control man-made vehicles [il foL [Tot llll [l2j]. 

Many types of mathematical models have been used 
to describe coherent swarms. One popular approach is 
based on a continuum approximation in which scalar and 
vector fields are used to describe all of the relevant quan- 
tities d, [IH, 14, 0, EH- Another popular approach is 
based on treating every biological or mechanical individ- 
ual as a discrete particle @, EI, 15, 17, 18j]. Depending 
on the problem, these individual-based models may be 
deterministic or stochastic. 

Regardless of the type of swarm model being used, one 
can see the emergence of ordered swarm states from an 
initial disordered state where individual particles have 
random velocity directions [13, 14, 17 1. These ordered 
states may be translational or rotational in motion, and 
they may be spatially distributed or localized in clusters. 

In particular, it is known that a localized swarm state 
may transition to a new dynamical region as the system 
parameters or the noise intensity is changed. For exam- 
ple, it has been shown in [l8| that a planar model of 
self-propelling particles interacting via a harmonic at- 
tractive potential in the presence of noise possesses a 
noise-induced transition whereby the translational mo- 
tion of the swarm breaks down into rotational motion. 

Another aspect of swarm modeling that has not yet 
been considered is the effect of time delayed interactions 
arising from finite communication times between indi- 
viduals. Much attention has been given to the effects of 
time delays in the context of physiology [l9| , optics [13] , 
neurons [2l|], lasers [23], and many other types of sys- 



tems. The aim of this Letter is to study the effect of a 
communication time delay on a model of self-propelling 
individuals that interact through a pairwise attractive 
force in the presence of noise. 

We consider a general two-dimensional (2D) model of 
a swarm that consists of identical self-propelled particles 
of unit mass. The model is described by the following 
evolution equations of motion: 
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where Ti(t) and Vi(t) are respectively the 2D position and 
velocity vectors of the z th particle at time t. The terms Vi 
and — |f i| 2 fi define respectively the mechanisms of self- 
propulsion and frictional drag. Therefore, if the last two 
terms on the right hand side of Eq. ([2]) are neglected, the 
particles will approach an equilibrium speed of v eq = 1 . 

The term Vi in Eq. ([2|) describes the social interaction, 
or communication, of the i th individual with all of the 
other individuals. There are many possible choices for 
Vi (e.g. Morse function, power law function, etc.). As an 
example, we define Vi as follows: 
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where a is the particle interaction coupling parameter, N 
is the number of particles, and r is a constant communi- 
cation time delay. This particular choice of V, assumes 
that only pairwise interactions are important. Further- 
more, the interaction is purely attractive and grows lin- 
early with the separation between two particles, much 
like a spring potential. 

Lastly, the term in Eq. @ describes a stochastic 
white force of intensity D. This noise is independent for 
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different particles, and is characterized by the following 

correlation functions: 

(&(*)) = 0, (&(t)fc(t)> = 2IW(t - f)V 

We numerically integrate Eqs. CO)-© using a stochas- 
tic fourth-order Runge-Kutta scheme with a constant 
time step size of 0.001. To achieve a traveling, local- 
ized swarm state, we used constant initial conditions Q 
and switched on the noise after a short amount of time 
had passed. 

It was shown in [18] that the model described by 
Eqs. O-Q with r = (i.e. no time delay) possesses 
a noise-induced transition whereby a large enough noise 
intensity causes a translating swarm of individuals to 
transition to a rotating swarm with a stationary cen- 
ter of mass, where the center of mass is defined as 
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FIG. 1: Snapshot of (a) a translating swarm (taken at t = 18), 
and (b) a rotating swarm (taken at t = 40), with a = 100, 
TV = 300, r = 0, and D = 0.08. 

Figures [lja) and QJb) show snapshots of a swarm at 
£ = 18 and t = 40 respectively, with a = 100, N = 300, 
t = 0, and D = 0.08. The noise was switched on at 
t = 10. One can see that the translating swarm [Fig.QJa)] 
has undergone a noise-induced breakdown to become a 
rotating swarm [Fig.QJb)]. For these values of a and N, 
if a noise intensity of D < 0.054 is used, then the swarm 
will continue to translate and it will not transition to a 



rotational state [24 1 . 

Regardless of which state the swarm is in (translating 
or rotating), the addition of a communication time de- 
lay leads to another type of transition. This transition 
occurs if the coupling parameter, a, is large enough. As 
an example, we consider a swarm that has already un- 
dergone a noise-induced transition to a rotational state 
before switching on the communication time delay. 

Figures [2ja)-[2](d) show snapshots of a swarm at t = 50, 
t = 100, t = 300, and t — 600 respectively, with a = 2, 
N — 300, t = 1, and D = 0.08. The noise was switched 
on at t — 10, and since the noise intensity, £>, is high 
enough, the noise caused the swarm to transition to a ro- 
tating state [similar to the one shown in Fig. [T][b)] . With 
the swarm in this stationary, rotating state, the commu- 
nication time delay was switched on at t — 40. One can 
see that for these values of time delay and coupling pa- 
rameter there is no qualitative change in the stationary, 
rotating swarm state. 
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FIG. 2: Snapshots of a swarm taken at (a) t — 50, (b) t — 100, 
(c) t = 300, and (d) t = 600, with a = 2, N = 300, and 
D — 0.08. The swarm was in a rotational state when the 
time delay of r = 1 was switched on at t — 40. 



In contrast to this, Figs. [3](a)-[3]j) show snapshots of 
a swarm at t = 50, t = 60, t — 62, t — 64, t = 66, 
t = 68, t = 70, t = 72, t = 74, and t = 76 respectively. 
As in the previous case, N — 300, r = 1, D = 0.08, 
the noise was switched on at t = 10 (causing the swarm 
to transition to a stationary, rotating state), and once 
in this rotating state, the time delay was switched on at 
t = 40. The only difference is that now the value of the 
coupling parameter is a — 4. One can see that with the 
evolution of time, the individual particles become aligned 
with one another and the swarm becomes more compact. 
Additionally, the swarm is no longer stationary, but has 
begun to oscillate [Figs. [3^g)-|3jj )] . This clockwise oscilla- 
tion can more clearly be seen in Fig. [H which consists of 
the center of mass, R, of the stationary, rotating swarm 
at t = 40 (denoted by a "cross" marker) along with snap- 
shots of the oscillating swarm taken at t = 90.2, t = 90.6, 
t = 91.0, and t = 91.4. 

This compact, oscillating aligned swarm state looks 
similar to a single "clump" that is described in [l(| • How- 
ever, where each "clump" of [l(| contains only some of 
the total number of swarming particles, our swarm con- 
tains every particle. Additionally, while a deterministic 
model along with global coupling is used to attain the 
"clumps" of [10( 1 , our oscillating aligned swarm is attained 
with the use of noise and a time delay. 

As we have shown, once the stochastic swarm is in the 
stationary, rotating state, the addition of a time delay in- 
duces an instability. At this point, the stochastic pertur- 
bations have a minimal effect on the swarm. Therefore, 
we will investigate the stability of the stationary, rotating 
swarm state by deriving the mean field equations with- 
out noise. The coordinates Xi and yi of each particle in 
the swarm can be written as follows: 



X + 5x t and ^ = Y + Syi 



(4) 
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FIG. 4: Motion of the oscillating swarm about the center of 
mass of the stationary, rotating swarm. The oscillating swarm 
is shown at t = 90.2 (left), t = 90.6 (top), t = 91.0 (right), 
and t = 91.4 (bottom). The location of the center of mass 
of the swarm (at t = 40) is denoted with a "cross" marker 
(center). 
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FIG. 3: Snapshots of a swarm taken at (a) t = 50, (b) t = 60, 
(c) t = 62, (d) t = 64, (e) t = 66, (f) t = 68, (g) t = 70, (h) 
t = 72, (i) t = 74, and (j) t = 76, with a = 4, N = 300, and 
D = 0.08. The swarm was in a rotational state when the time 
delay of r = 1 was switched on at t — 40. 



where X and Y are the coordinates of the center of mass, 
R, of the swarm. Substitution of Eq. ((4|) into the second- 
order differential equation that is equivalent to Eqs. (fT|)- 
(|3"|) gives an evolution equation for each Xi and i/j. Sum- 
ming all i of these equations, using the fact that 

N N 

-^2xi{t)=X(t) and -J2vi(t)=Y(t), (5) 



and ignoring all fluctuation terms, leads to the following 
zero-order mean field equations for the center of mass: 



X(t) = [{l-X 2 )X-XY 2 }(t)-a{X{t)-X(t-T)) : (6) 



Y(t) = [(1 - Y 2 )Y - YX 2 ]{t) - a(Y(t) - Y(t - r)). (7) 

The steady state is given by X(t) = Y(t) = 0, 
X(t) = X(t - t), and Y(t) = Y(t - r). Consideration of 
small disturbances about the steady state allows one to 
determine the linear stability. The characteristic equa- 
tion associated with the linearization of Eqs. ©-Q is 



A(l-A)+ae 



-At 



a = 0, 



(8) 



where the exponential term exp(— At) is due to the time 
delay in the governing equations. Since Eq. ([8]) is tran- 
scendental (which is often the case for delay differential 
equations), there exists the possibility of an infinite num- 
ber of solutions. 

Our numerical simulations indicate the existence of a 
supercritical Hopf bifurcation as the value of the cou- 
pling parameter, a, is increased (Figs. EE])- We identify 
the Hopf bifurcation point by choosing the eigenvalue to 
be purely imaginary. Then our choice of A = iu> is substi- 
tuted into Eq. (jHJl . The separation of Eq. ((SJ) into real and 
imaginary parts leads to an equation for the frequency, 
u>, along with an equation for the value of a at the Hopf 
bifurcation point. The two equations are 



uj 2 + lo cot (wt) — lo csc (wt) = 0, 



sin (lot) 



(9) 
(10) 



Given a specific value of r, Eq. |9]) can be solved nu- 
merically for uj. These values of r and u> can then be 
substituted into Eq. (TTTj| to determine the value of a at 
the Hopf point. 



4 



so 



«P= 40 



15.2 


















15 


y>v i=6s 






/ t=0.25 


- ^14.8 










14.6 
14.4 




/ /J/ 1=72 


/ T=0.3 






15.1 15.3 15.5 


15.7 15.9 












/ T=0.35 








T=0.5 






T=1.5 











2 4 6 8 10 12 14 



FIG. 5: Comparison of analytical (solid line) and numerical 
("cross" markers) values of an and w for several choices of r. 
The analytical result is found using Eqs. ©-([TO]), while the 
numerical result is found using a continuation method [25| for 
Eqs. (©-(0. The inset shows the stochastic trajectory of the 
center of mass of the swarm from t = 45 to t — 90 for the 
example shown in Fig. [3] 

Figure [5] shows an excellent comparison of the analyt- 
ical result given by Eqs. (|9|)- ([T0|) with a numerical result 
which was found using a continuation method [25I ] for 
Eqs. ©-(O for several choices of r. Furthermore, for 
t = 1, the value of a at the bifurcation point is an ~ 3.2. 
This value of an corresponds very well to the change 
in behavior of the stochastic swarm that was seen as the 
value of the coupling parameter was increased from a = 2 
to a = 4 (Figs. HP . 

More evidence of the Hopf bifurcation is seen in the 
inset of Fig. [51 The inset shows the stochastic trajectory 
of the center of mass of the swarm from t = 45 to t = 90 
for the example shown in Fig. [3l Once the time delay 
is switched on at t — 40 (with the swarm located at the 
center of the inset figure), the swarm begins to oscillate. 
The swarm moves along an elliptical path [the position 
of its center of mass is denoted at several times that cor- 
respond to Figs. G2», 13(d), Hf),Hh), and Eft)], until it 
eventually converges to the circular limit cycle. 

To summarize, we studied the dynamics of a self- 
propelling swarm in the presence of noise and a constant 
communication time delay and prove that the delay in- 
duces a transition that depends upon the size of the inter- 
action coupling coefficient. Although our analytical and 
numerical results were obtained using a model with lin- 
ear, attractive interactions, the analysis may be applied 
to models with more general forms of social interaction 
(these results will appear elsewhere). 

Our results provide insight into the stability of com- 
plex systems comprised of individuals interacting with 
one another with a finite time delay in a noisy environ- 
ment. Furthermore, the results may prove to be useful 
in controlling man-made vehicles where actuation and 



communication are delayed, as well as in understanding 
swarm alignment in biological systems. 
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